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LONG-TERM  GOALS 

Smoothed  Particle  Hydrodynamics  (SPH)  is  a  meshless  numerical  method  that  is  being  developed  for 
the  study  of  nearshore  waves  and  other  Navy  needs.  The  Lagrangian  nature  of  SPH  allows  the 
modeling  of  wave  breaking,  surf  zones,  ship  waves,  and  wave-structure  interaction,  where  the  free 
surface  becomes  convoluted  or  splash  occurs. 

OBJECTIVES 

To  improve  the  ability  of  the  meshfree  Lagrangian  numerical  method  Smoothed  Particle  Hydro¬ 
dynamics  (SPH)  to  be  a  useful  hydrodynamics  model  for  breaking  waves  and  the  nearshore  zone, 
particularly  for  case  where  spray  and  splash  are  important.  To  utilize  the  massively  parallel  graphics 
processing  units  on  computers  to  develop  the  GPU-accelerated  model  GPUSPH  to  solve  a  number  of 
problems  relevant  to  the  U.S.  Navy.  The  science  objective  is  to  be  able  to  accurately  model  the 
complex  flows  associated  with  breaking  water  waves,  including  instantaneous  motions  as  well  as 
(time-averaged)  wave-induced  flows,  such  as  undertow,  longshore  currents,  and  rip  currents. 

APPROACH 

The  approach  is  based  on  improving  various  aspects  of  the  SPH  code,  including  the  development  of  a 
multi-graphics  processing  unit  (GPU)  version  of  the  code  (GPUSPH);  applying  the  code  to  more 
validation  tests;  and  to  examine  in  some  detail  new  aspects  of  the  model  by  applying  it  to  different 
situations  relevant  to  Navy  needs. 

WORK  COMPLETED 

•  Open  source  version  of  GPUSPH  released  on  the  web:  January  17,  201 1 

•  GPUSPH  model  was  extended  to  find  fluid  parameters  at  fixed  (Eulerian)  test  points 

•  Further  validation  of  nearshore  modeling,  including  waves,  a  rip  current  system  and  related  surf 
zone  circulation  using  the  test  points 
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Free-surface  particle  detection  implemented  to  simplify  determination  of  the  free  surface  location 


RESULTS 

Open  source  release— the  GPUSPH  model  was  released  as  open  source  in  January  after  a  month  long 
interaction  with  some  of  the  code  developers  (Alexis  Herault,  Eugenio  Rustico)  in  Baltimore, 
Maryland.  The  code,  in  C++,  with  test  problems,  is  available  at  www.ce.jhu.edu/dalrymple/GPUSPH. 

Test  particles— The  Lagrangian  nature  of  SPH  acts  as  an  advantage  to  model  problems  with  large 
deformations.  On  the  other  hand,  as  the  SPH  particles  (or  nodes)  are  moving  with  each  time  step,  it  is 
difficult  to  record  the  fluid  parameters  in  fixed  positions  as  might  be  recorded  in  a  field  or  laboratory 
experiment.  In  order  to  look  at  time-averaged  horizontal  velocity  profiles  and  time-averaged  water 
elevations,  velocities  and  elevations  are  required  to  be  measured  at  fixed  spatial  positions. 

To  solve  this  problem,  the  GPUSPH  was  extended  to  find  fluid  parameters  at  fixed  spatial  positions  (as 
well  on  the  moving  particles).  A  set  of  testpoints,  which  have  fixed  positions  in  the  computational 
domain,  were  introduced.  Velocity  and  pressure  at  these  testpoints  are  computed  as  a  function  of  time 
with  the  SPH  interpolation  formula,  using  the  properties  of  moving  particles  that  are  currently  located 
in  the  neighborhood  of  the  testpoint. 

In  order  to  validate  the  extended  model,  experimental  results  of  a  rip  current  system  (Dronen  et  al., 
2002)  have  been  compared  to  the  numerical  results  of  extended  GPUSPH  model.  The  laboratory 
experiments  were  done  in  a  4  m  wide  and  30  m  long  wave  tank.  An  idealized  sand  bar  was  located  on 
one  side  of  the  tank  that  triggers  local  breaking,  leading  to  the  rip  current  creation.  Figure  (1)  shows  the 
numerical  results.  Rip  currents,  which  are  strong  currents  flowing  in  the  oposite  direction  of  water 
waves,  are  evident  on  the  right  side  of  wave  tank. 
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Figure  1.  Numerical  results  of  rip  current  simulation  based  on  the  laboratory  experiments  of 
Dronen  et  al.  2001.  The  wavemaker  is  located  at  the  upper  right,  and  an  idealized  sand  bar  on  the 
left  side.  In  this  frame,  wave  breaking,  harmonic  generation  over  the  bar,  wave  refraction, 
diffraction,  shoaling,  and  wave-current  interaction  are  all  visible.  Visualization  provided  by 

Templeman  Automation. 

The  numerical  model  is  capable  of  modelling  wave-current  interaction,  harmonic  generation  (a  shallow 
water  nonlinear  process),  wave  breaking,  wave  diffraction,  refraction  and  nearshore  circulations.  The 
rip  current  system,  which  consists  of  (1)  an  offshore-directed  rip  current  in  the  channel  on  the  right  of 
Figure  1,  (2)  a  circulation  cell  over  the  sand  bar,  and  (3)  an  onshore  oppositely-rotating  circulation  cell 
in  the  bar  trough  is  sucessfully  modeled. 

Figure  2  shows  the  mean  vorticity,  mean  horizontal  velocity  profile  and  related  wave-induced 
circulation.  Here,  wavemaker  is  located  outside  of  the  figure  to  the  right  and  beach  is  located  in  the  left 
hand  side.  In  the  bottom  frame  of  Figure  2,  the  mean  water  elevation  is  also  shown  in  different  parts  of 
the  tank.  A  large  portion  of  the  elevation  is  the  set-up  on  the  beach  and  in  the  bar  trough,  incuded  by 
the  breaking  on  the  bar.  Gradients  of  the  mean  water  level  provide  the  hydrostatic  head  to  drive  the  rip 
current  (in  the  rip  channel)  and  the  associated  circulation  cells. 
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Figure  2.  Top:  Mean  vorticity  associated  with  the  wave-induced  currents  and  mean  horizontal 
velocity.  Bottom:  Mean  water  level  elevations,  including  the  wave  set-up. 

GPU-SPH  is  capable  of  three  dimentional  analysis  of  flow  which  leads  to  a  better  understanding  of 
flow  structure.  Three  dimentional  velocity  profiles  in  different  positions  located  in  rip  channel,  over  the 
bar  and  near  the  beach,  are  illustrated  in  figure  (3). 
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Figure  3.  Three  dimensional  velocity  profiles  at  fixed  locations  over  the  Dronen  et  al.  bathymetry 


Free-surface  particle  detection— To  study  free  surface  flows  and  analyze  their  complex  deformations, 
we  need  to  know  which  particles  are  located  on  the  free  surface.  Free  surface  detection  algorithm  that 
is  used  here  is  similar  to  Marrone  et  al.  (2010)  approach  and  consists  of  two  steps:  first,  the  local 
normal  vector  for  each  particle  is  computed  using  gradient  of  kernel  function.  In  the  second  step, 
geometric  properties  of  free  surface  particles  are  used  for  their  accurate  detection.  For  each  particle  a 
cone  is  defined  with  particle’s  normal  vector  as  its  axis  and  a  cone  angle  that  has  been  set  to  be  equal  to 
.  Then  a  control  check  is  applied  to  find  whether  or  not  at  least  one  neighboring  particle  exist  in  this 
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cone  region.  If  no  neighboring  particle  is  found  within  the  cone  region,  then  that  particle  belongs  to  the 
set  of  surface  particles.  This  control  condition  is  given  as: 

/7.F. 

V j, — —<  cos  {cone  _angle)  =>  /  E  Free  _  surface  _  particles 
r 

where  i  is  the  particle  of  interest,  j  is  its  neighboring  particle,  n  is  local  normal  vector  for  particle  i 
and  r  is  the  distance  from  particle  i  to  particle  j .  Figure  (4)  shows  a  sketch  of  fluid  domain  and  free 

surface  detection  using  cone  region.  Since  we  are  using  the  gradient  of  kernel  function  that  has  been 
computed  previously  in  another  part  of  the  model,  this  algorithm  doesn’t  impose  a  significant  increase 
in  computational  cost. 
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Figure  4.  Schematic  plan  of fluid  domain  and  free  surface 


The  free  surface  detection  algorithm  as  applied  to  different  cases  is  shown  in  figure  5  and  6.  Free 
surface  particles  are  shown  in  red  and  other  particles  are  shown  in  light  blue.  Note  that  the  model 
resolution  is  deliberately  low,  so  that  the  free  surface  particles  are  readily  viewed.  The  first  case  shows 
a  plunging  wave.  As  it  can  be  seen  in  the  picture,  free  surface  detection  algorithm  is  capable  of 
capturing  the  cavity  that  occurs  with  a  plunging  wave.  The  second  case  shows  a  dam  break  problem 
and  the  interaction  of  the  waves  with  a  structure.  The  last  case  shows  an  open  channel  problem  with 
periodic  boundary  condition. 
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Figure  5.  Free  surface  detection  for  plunging  wave.  The  irregularity  of  the  free  surface  particles  is 
due  to  the  low  number  of particles  used  to  clearly  show  the  free  surface  particles,  and  the  3D  nature 

of  the  irregular  free  surface  as  viewed  from  the  side. 
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Figure  6.  Free  surface  detection  for  (a)  dam  break  problem  with  an  obstacle, 
and  (b)  open  channel  flow  problem 

International  Collaborations:  We  currently  have  a  strong  international  collaboration  in  the  developer 
group  for  GPUSPH  with  the  Istituto  Nazionale  di  Geofisica  e  Vulcanologia  (sezione  di  Catania),  the 
Universita  di  Catania,  Conservatoire  National  des  Arts  et  Metiers,  Paris,  for  the  development  of  GPU- 
SPHysics.  Dr.  Herault  (CNAM  and  INGV),  Dr.  Bilotta  (UC),  and  Mr.  Eugenio  Rustico  (UC)  are 
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members  of  this  collaboration.  A  new  collaborative  organization  is  being  set  up  in  October  2011  that 
will  include  the  original  GPUSPH  developers  and  EDF  R&D  (the  research  and  development  division 
of  Electricite  de  France)  for  the  joint  development  of  GPUSPH.  EDF  has  an  existing  SPH  code, 
SPARTACUS,  but  is  interested  in  co-developing  the  GPU  accelerated  GPUSPH  code.  It  is  anticipated 
that  other  European  institutes  will  join  the  effort. 

Ongoing  efforts — There  are  three  major  initiatives  underway.  The  first  involves  adding  the  ability  to 
model  floating  bodies  with  GPUSPH.  This  involves  determining  the  forces  on  an  object  and  then 
allowing  the  object  to  respond  to  the  forces  appropriately.  The  second  is  lead  by  Eugenio  Rustico, 
which  is  the  development  of  multi-GPU  configurations  of  the  GPUSPH.  This  is  critical  as  the 
GPUSPH  code  on  a  single  GPU  card  is  only  capable  of  handling  5  million  particles  in  a  simulation. 
JHU  is  currently  completing  the  construction  of  a  100  GPU  computing  cluster,  which  will  permit 
extensive  testing  of  high  resolution  GPUSPH  simulations.  Finally  multi-fluid  versions  of  GPUSPH 
are  underway. 

IMPACT/APPLICATIONS 

Smoothed  Particle  Hydrodynamics  is  proving  to  be  a  competent  modeling  scheme  for  free  surface 
flows  in  two  and  three  dimensions.  As  the  GPU  hardware  improves,  it  is  expected  that  the  resolution 
of  SPH  will  increase  tremendously  bringing  the  modeling  into  realistically  sized  domains. 
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